Exercício de Métodos Inferenciais Avançados 2/18
Há 41 observações com alguma informação faltando e 1409 casos completos.
alias(lm(BirthWeightOz ~ . - BirthWeightGm, data = new_ncbirths(NCbirths)))## Model :
## BirthWeightOz ~ (Plural + Sex + MomAge + Weeks + Marital + RaceMom +
## HispMom + Gained + Smoke + BirthWeightGm + Low + Premie +
## MomRace) - BirthWeightGm
##
## Complete :
## (Intercept) PluralTwins PluralTriplets SexF MomAge Weeks
## MomRacehispanic 0 0 0 0 0 0
## MomRaceother 0 0 0 0 0 0
## MomRacewhite 1 0 0 0 0 0
## MaritalNot Married RaceMomBlack RaceMomAmericanIndian
## MomRacehispanic 0 0 0
## MomRaceother 0 0 1
## MomRacewhite 0 -1 -1
## RaceMomChinese RaceMomHispanic RaceMomFilipino
## MomRacehispanic 0 1 0
## MomRaceother 1 0 1
## MomRacewhite -1 -1 -1
## RaceMomOtherAsianOrPacific HispMomM HispMomN HispMomO
## MomRacehispanic 0 0 0 0
## MomRaceother 1 0 0 0
## MomRacewhite -1 0 0 0
## HispMomP HispMomS Gained SmokeY LowY PremieY
## MomRacehispanic 0 0 0 0 0 0
## MomRaceother 0 0 0 0 0 0
## MomRacewhite 0 0 0 0 0 0
alias(lm(BirthWeightOz ~ . - BirthWeightGm, data = new_ncbirths(NCbirths, T)))## Model :
## BirthWeightOz ~ (Plural + Sex + MomAge + Weeks + Marital + Gained +
## Smoke + BirthWeightGm + Low + Premie + Etnicidade) - BirthWeightGm
A sobreposição de informações dadas nas colunas étnicas aparece no relatório criado pela função alias(). Isso pôde ser corrigido ao se unificar as informações em uma única coluna.
Curiosamente, há combinações inesperadas ali, como mães negras oriundas de Porto Rico e da América do Sul sendo consideradas também hispânicas. Por outro lado, distinguir índios de continentes distintos parece compreensível.
Da mesma forma que no modelo 1, vejamos com o que se parece um modelo com “tudo” dentro. Entretanto, como há campos com forte correlação entre si e com a variável dependente, precisamos nos livrar antes de:
Além disso, os 40 casos de ganho de peso ausentes precisam ser avaliados para uma tomada de decisão.
# ncb_weight_noNAs = complete.cases(ncb_weight)
# fit_weight_1 é o fModelo 3 ====
fit_weight_1 = lm(BirthWeightOz ~ ., data=ncb_weight, na.action = na.omit)##
## Call:
## lm(formula = BirthWeightOz ~ ., data = ncb_weight, na.action = na.omit)
##
## Residuals:
## Min 1Q Median 3Q Max
## -60.588 -10.165 -0.473 10.194 51.977
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -31.76600 10.41147 -3.051 0.002324 **
## PluralTwins -24.91791 2.70283 -9.219 < 2e-16 ***
## PluralTriplets -32.97522 8.37525 -3.937 8.65e-05 ***
## SexF -3.29641 0.87689 -3.759 0.000178 ***
## MomAge 0.38212 0.08229 4.644 3.75e-06 ***
## Weeks 3.46722 0.24046 14.419 < 2e-16 ***
## MaritalNot Married -1.84540 1.14051 -1.618 0.105881
## Gained 0.28078 0.03227 8.701 < 2e-16 ***
## SmokeY -7.09110 1.29607 -5.471 5.30e-08 ***
## PremieY -7.48998 1.88547 -3.972 7.48e-05 ***
## Etnicidadeblack Black N -2.35540 3.70339 -0.636 0.524873
## Etnicidadeblack portoriq 7.33369 16.86260 0.435 0.663696
## Etnicidadeblack south american 4.18383 16.81073 0.249 0.803492
## Etnicidadecentro-south american 0.49574 5.02110 0.099 0.921366
## Etnicidadechinese 3.17904 12.16810 0.261 0.793931
## Etnicidadecuban 3.71059 12.14769 0.305 0.760064
## Etnicidadefilipino -30.03061 16.81461 -1.786 0.074320 .
## Etnicidadehispanic other 3.68649 10.15569 0.363 0.716662
## Etnicidademexican 2.36859 3.91357 0.605 0.545129
## EtnicidadeOtherAsianOrPacific -1.72621 5.03311 -0.343 0.731673
## Etnicidadeportoriq -13.26520 7.18070 -1.847 0.064911 .
## Etnicidadesouth americanIndian 7.88998 16.81736 0.469 0.639031
## Etnicidadewhite 1.41976 3.65717 0.388 0.697918
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.39 on 1386 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.4597, Adjusted R-squared: 0.4512
## F-statistic: 53.61 on 22 and 1386 DF, p-value: < 2.2e-16
##
## Call:
## lm(formula = BirthWeightOz ~ ., data = ncb_weigh_patchedNA, na.action = na.omit)
##
## Residuals:
## Min 1Q Median 3Q Max
## -60.250 -10.381 -0.569 10.246 52.373
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -36.91916 10.24067 -3.605 0.000323 ***
## PluralTwins -23.90895 2.65762 -8.996 < 2e-16 ***
## PluralTriplets -32.33524 8.42338 -3.839 0.000129 ***
## SexF -3.41657 0.87099 -3.923 9.17e-05 ***
## MomAge 0.35694 0.08148 4.381 1.27e-05 ***
## Weeks 3.62095 0.23489 15.415 < 2e-16 ***
## MaritalNot Married -2.40859 1.12750 -2.136 0.032830 *
## Gained 0.28113 0.03244 8.666 < 2e-16 ***
## SmokeY -6.77780 1.29182 -5.247 1.78e-07 ***
## PremieY -6.62059 1.86564 -3.549 0.000400 ***
## Etnicidadeblack Black N -2.34155 3.72682 -0.628 0.529910
## Etnicidadeblack portoriq 7.00896 16.97642 0.413 0.679768
## Etnicidadeblack south american 4.47865 16.92608 0.265 0.791355
## Etnicidadecentro-south american 1.13483 5.00139 0.227 0.820532
## Etnicidadechinese 3.31510 12.25095 0.271 0.786738
## Etnicidadecuban 3.89491 12.23136 0.318 0.750201
## Etnicidadefilipino -29.98661 16.92982 -1.771 0.076736 .
## Etnicidadehispanic other 3.80635 10.22475 0.372 0.709748
## Etnicidademexican 2.49703 3.91547 0.638 0.523749
## EtnicidadeOtherAsianOrPacific -2.40474 5.01479 -0.480 0.631635
## Etnicidadeportoriq -9.51223 6.88705 -1.381 0.167441
## Etnicidadesouth americanIndian 7.62941 16.93276 0.451 0.652367
## Etnicidadewhite 1.11467 3.68020 0.303 0.762023
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.5 on 1426 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.4619, Adjusted R-squared: 0.4536
## F-statistic: 55.65 on 22 and 1426 DF, p-value: < 2.2e-16
As diferenças entre os coeficientes não são estatisticamente significantes:
comparaCoeficientes(sf1, sf1.1, "Coeficientes com e sem linhas incompletas")Chama a atenção o coeficiente para Etnicidadefilipino, com valor absoluto elevado e com significância de 0.055637 0.074320 (0.076736 no modelo que inclui as linhas com valores ausentes), perto do limite escolhido para \(\alpha\). Existe exatamente 1 registro de mãe com essa etnia.
if (useEtn) {
ncb_weight_NF = ncb_weight[ncb_weight$Etnicidade != "filipino",]
} else {
ncb_weight_NF = ncb_weight[ncb_weight$RaceMom != "Filipino",]
}
fit_weight_1b = lm(BirthWeightOz ~ ., data=ncb_weight_NF, na.action = na.omit)
sf1b = summary(fit_weight_1b)
sf1b##
## Call:
## lm(formula = BirthWeightOz ~ ., data = ncb_weight_NF, na.action = na.omit)
##
## Residuals:
## Min 1Q Median 3Q Max
## -60.588 -10.192 -0.478 10.198 51.977
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -31.76600 10.41147 -3.051 0.002324 **
## PluralTwins -24.91791 2.70283 -9.219 < 2e-16 ***
## PluralTriplets -32.97522 8.37525 -3.937 8.65e-05 ***
## SexF -3.29641 0.87689 -3.759 0.000178 ***
## MomAge 0.38212 0.08229 4.644 3.75e-06 ***
## Weeks 3.46722 0.24046 14.419 < 2e-16 ***
## MaritalNot Married -1.84540 1.14051 -1.618 0.105881
## Gained 0.28078 0.03227 8.701 < 2e-16 ***
## SmokeY -7.09110 1.29607 -5.471 5.30e-08 ***
## PremieY -7.48998 1.88547 -3.972 7.48e-05 ***
## Etnicidadeblack Black N -2.35540 3.70339 -0.636 0.524873
## Etnicidadeblack portoriq 7.33369 16.86260 0.435 0.663696
## Etnicidadeblack south american 4.18383 16.81073 0.249 0.803492
## Etnicidadecentro-south american 0.49574 5.02110 0.099 0.921366
## Etnicidadechinese 3.17904 12.16810 0.261 0.793931
## Etnicidadecuban 3.71059 12.14769 0.305 0.760064
## Etnicidadehispanic other 3.68649 10.15569 0.363 0.716662
## Etnicidademexican 2.36859 3.91357 0.605 0.545129
## EtnicidadeOtherAsianOrPacific -1.72621 5.03311 -0.343 0.731673
## Etnicidadeportoriq -13.26520 7.18070 -1.847 0.064911 .
## Etnicidadesouth americanIndian 7.88998 16.81736 0.469 0.639031
## Etnicidadewhite 1.41976 3.65717 0.388 0.697918
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.39 on 1386 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.4594, Adjusted R-squared: 0.4513
## F-statistic: 56.1 on 21 and 1386 DF, p-value: < 2.2e-16
comparaCoeficientes(sf1b, sf1)vif(fit_weight_1)## GVIF Df GVIF^(1/(2*Df))
## Plural 1.173561 2 1.040822
## Sex 1.007636 1 1.003811
## MomAge 1.315787 1 1.147078
## Weeks 2.134540 1 1.461006
## Marital 1.541343 1 1.241509
## Gained 1.050612 1 1.024994
## Smoke 1.099646 1 1.048640
## Premie 2.087255 1 1.444734
## Etnicidade 1.419620 13 1.013568
vif(fit_weight_1b)## GVIF Df GVIF^(1/(2*Df))
## Plural 1.173533 2 1.040816
## Sex 1.006888 1 1.003438
## MomAge 1.313160 1 1.145932
## Weeks 2.134147 1 1.460872
## Marital 1.540767 1 1.241276
## Gained 1.050590 1 1.024983
## Smoke 1.099513 1 1.048576
## Premie 2.087037 1 1.444658
## Etnicidade 1.415500 12 1.014584
As colunas Weeks e Premie são as que mais inflacionam a variância desse modelo. A retirada da mãe filipina diminuiu levemente a inflação do campo Etnicidade.
Menos é mais (de novo)
Baseado no contexto e nos valores observados das significâncias dos coeficientes do modelo 3, seleciono um conjunto menor de colunas para incluir no modelo.
#COLSET_01 = c("Plural", "Sex", "MomAge", "Weeks",
# "Gained", "Smoke",
# "BirthWeightOz")
# fit_weight_2 é o fModelo 4 ====
fit_weight_2 = lm(
BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke,
data = ncb_weight
)
summary(fit_weight_2)##
## Call:
## lm(formula = BirthWeightOz ~ Plural + Sex + MomAge + Weeks +
## Gained + Smoke, data = ncb_weight)
##
## Residuals:
## Min 1Q Median 3Q Max
## -63.672 -10.236 -0.157 10.486 49.522
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -63.14127 7.06144 -8.942 < 2e-16 ***
## PluralTwins -25.54412 2.71720 -9.401 < 2e-16 ***
## PluralTriplets -33.92719 8.45212 -4.014 6.28e-05 ***
## SexF -3.42455 0.88540 -3.868 0.000115 ***
## MomAge 0.50662 0.07338 6.904 7.65e-12 ***
## Weeks 4.16276 0.17797 23.390 < 2e-16 ***
## Gained 0.28536 0.03215 8.875 < 2e-16 ***
## SmokeY -7.24103 1.26076 -5.743 1.14e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.59 on 1401 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.4408, Adjusted R-squared: 0.438
## F-statistic: 157.8 on 7 and 1401 DF, p-value: < 2.2e-16
fit_weight_2b = lm(
BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke,
data = ncb_weight_NF)
summary(fit_weight_2b)##
## Call:
## lm(formula = BirthWeightOz ~ Plural + Sex + MomAge + Weeks +
## Gained + Smoke, data = ncb_weight_NF)
##
## Residuals:
## Min 1Q Median 3Q Max
## -63.654 -10.185 -0.172 10.478 49.537
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -63.40360 7.05644 -8.985 < 2e-16 ***
## PluralTwins -25.56954 2.71478 -9.419 < 2e-16 ***
## PluralTriplets -33.98035 8.44453 -4.024 6.03e-05 ***
## SexF -3.38139 0.88489 -3.821 0.000139 ***
## MomAge 0.51259 0.07339 6.985 4.39e-12 ***
## Weeks 4.16531 0.17782 23.425 < 2e-16 ***
## Gained 0.28563 0.03213 8.891 < 2e-16 ***
## SmokeY -7.25583 1.25964 -5.760 1.03e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.57 on 1400 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.4419, Adjusted R-squared: 0.4391
## F-statistic: 158.3 on 7 and 1400 DF, p-value: < 2.2e-16
anova(fit_weight_1, fit_weight_2)## Analysis of Variance Table
##
## Model 1: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Marital + Gained +
## Smoke + Premie + Etnicidade
## Model 2: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1386 372384
## 2 1401 385454 -15 -13069 3.2429 2.512e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(fit_weight_1b, fit_weight_2b)## Analysis of Variance Table
##
## Model 1: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Marital + Gained +
## Smoke + Premie + Etnicidade
## Model 2: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1386 372384
## 2 1400 384483 -14 -12098 3.2164 4.908e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Inesperadamente, o a soma dos quadrados dos resíduos aumentou quando as outras variáveis foram retiradas.
# fit_weight_5 é o fModelo 5 ====
if (useEtn) {
fit_weight_5 = lm(
BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke + Etnicidade,
data=ncb_weight)
sf5 = summary(fit_weight_5)
print(sf5)
anova(fit_weight_1, fit_weight_2, fit_weight_5)
} else {
fit_weight_5 = lm(
BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke + RaceMom,
data = ncb_weight)
sf5 = summary(fit_weight_5)
print(sf5)
# fit_weight_5{b,c} são variações do fModelo 3 ====
fit_weight_5b = lm(
BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke + MomRace,
data = ncb_weight)
summary(fit_weight_5b)
fit_weight_5c = lm(
BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke + HispMom,
data = ncb_weight)
summary(fit_weight_5c)
anova(fit_weight_1, fit_weight_2, fit_weight_5, fit_weight_5b, fit_weight_5c)
}##
## Call:
## lm(formula = BirthWeightOz ~ Plural + Sex + MomAge + Weeks +
## Gained + Smoke + Etnicidade, data = ncb_weight)
##
## Residuals:
## Min 1Q Median 3Q Max
## -65.805 -10.247 -0.376 10.312 52.716
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -61.07108 7.88766 -7.743 1.87e-14 ***
## PluralTwins -26.01446 2.70672 -9.611 < 2e-16 ***
## PluralTriplets -35.03222 8.40840 -4.166 3.29e-05 ***
## SexF -3.33780 0.88213 -3.784 0.000161 ***
## MomAge 0.44590 0.07632 5.842 6.40e-09 ***
## Weeks 4.12488 0.17780 23.199 < 2e-16 ***
## Gained 0.27935 0.03246 8.606 < 2e-16 ***
## SmokeY -7.66939 1.28020 -5.991 2.66e-09 ***
## Etnicidadeblack Black N -2.18985 3.72465 -0.588 0.556672
## Etnicidadeblack portoriq 11.56560 16.92197 0.683 0.494426
## Etnicidadeblack south american 5.59946 16.90440 0.331 0.740511
## Etnicidadecentro-south american 0.66689 5.05104 0.132 0.894979
## Etnicidadechinese 5.60070 12.22582 0.458 0.646949
## Etnicidadecuban 4.74089 12.21830 0.388 0.698064
## Etnicidadefilipino -29.30335 16.91408 -1.732 0.083410 .
## Etnicidadehispanic other 5.84796 10.20114 0.573 0.566558
## Etnicidademexican 2.61558 3.93672 0.664 0.506541
## EtnicidadeOtherAsianOrPacific -1.46203 5.05622 -0.289 0.772506
## Etnicidadeportoriq -12.73123 7.20327 -1.767 0.077377 .
## Etnicidadesouth americanIndian 9.05507 16.91154 0.535 0.592433
## Etnicidadewhite 2.50894 3.65806 0.686 0.492912
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.49 on 1388 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.4524, Adjusted R-squared: 0.4445
## F-statistic: 57.33 on 20 and 1388 DF, p-value: < 2.2e-16
## Analysis of Variance Table
##
## Model 1: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Marital + Gained +
## Smoke + Premie + Etnicidade
## Model 2: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke
## Model 3: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke +
## Etnicidade
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1386 372384
## 2 1401 385454 -15 -13069.3 3.2429 2.512e-05 ***
## 3 1388 377445 13 8008.8 2.2930 0.00536 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
comparaCoeficientes(summary(fit_weight_5), sf1)A inclusão dos campos étnicos prejudicou o modelo, pois os erros cresceram ou ficaram confusos (p>0.05).
Houve diferença significativa entre os coeficientes do intercepto e do campo Weeks (mas não sei por que)
O modelo 5 chamou atenção para o campo Weeks. Pelo contexto, espera-se uma forte correlação entre o peso do nascituro e a idade gestacional do parto: quanto mais próximo do termo, maior o peso (o que acontece depois do termo?).
Esta seção tenta isolar a influência desse campo no modelo.
# fit_weight_W é um modelo linear simples para comparar com o fModelo 6 ====
fit_weight_W = lm(BirthWeightOz ~ Weeks, data=ncb_weight)
sfW = summary(fit_weight_W)
print(sfW)##
## Call:
## lm(formula = BirthWeightOz ~ Weeks, data = ncb_weight)
##
## Residuals:
## Min 1Q Median 3Q Max
## -68.480 -11.994 -0.286 11.908 58.006
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -73.1171 6.7817 -10.78 <2e-16 ***
## Weeks 4.9028 0.1752 27.99 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 17.99 on 1447 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.3512, Adjusted R-squared: 0.3508
## F-statistic: 783.4 on 1 and 1447 DF, p-value: < 2.2e-16
Um modelo que tem apenas o campo Weeks como preditor de BirthWeighOz explica pouco mais de 35% da variância. É bastante, comparado com o modelo 3, com todas as variáveis, que explica 45.9741839%
# anova(fit_weight_1, fit_weight_W, fit_weight_2)
print("Cadê o noNA?")## [1] "Cadê o noNA?"
if (useEtn) { # Campo "etnicidade" unificado no modelo 6 ====
fit_weight_6 = lm(BirthWeightOz ~ Plural + Sex + MomAge + Gained + Smoke + Etnicidade, data=ncb_weight)
sf6 = summary(fit_weight_6)
print(sf6)
anova(fit_weight_1, fit_weight_6) # compara variância de "1" com "6"
} else { # Campos étnicos separados (como no original) no modelo 6 ====
fit_weight_6 = lm(BirthWeightOz ~ Plural + Sex + MomAge + Gained + Smoke + RaceMom, data = ncb_weight)
sf6 = summary(fit_weight_6)
print(sf6)
fit_weight_6b = lm(BirthWeightOz ~ Plural + Sex + MomAge + Gained + Smoke + MomRace, data = ncb_weight)
sf6b = summary(fit_weight_6b)
print(sf6b)
fit_weight_6c = lm(BirthWeightOz ~ Plural + Sex + MomAge + Gained + Smoke + HispMom, data = ncb_weight)
sf6c = summary(fit_weight_6c)
print(sf6c)
anova(fit_weight_1, fit_weight_2, fit_weight_6, fit_weight_6b, fit_weight_6c)
}##
## Call:
## lm(formula = BirthWeightOz ~ Plural + Sex + MomAge + Gained +
## Smoke + Etnicidade, data = ncb_weight)
##
## Residuals:
## Min 1Q Median 3Q Max
## -94.525 -10.042 0.959 11.864 52.319
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 92.92254 5.01749 18.520 < 2e-16 ***
## PluralTwins -45.55676 3.02917 -15.039 < 2e-16 ***
## PluralTriplets -67.39100 9.76461 -6.902 7.80e-12 ***
## SexF -2.27666 1.03741 -2.195 0.0284 *
## MomAge 0.54296 0.08974 6.050 1.86e-09 ***
## Gained 0.35434 0.03804 9.316 < 2e-16 ***
## SmokeY -8.92184 1.50624 -5.923 3.97e-09 ***
## Etnicidadeblack Black N -2.51296 4.38614 -0.573 0.5668
## Etnicidadeblack portoriq 8.28203 19.92673 0.416 0.6777
## Etnicidadeblack south american 5.05169 19.90672 0.254 0.7997
## Etnicidadecentro-south american -0.57218 5.94780 -0.096 0.9234
## Etnicidadechinese 1.43880 14.39566 0.100 0.9204
## Etnicidadecuban 6.84325 14.38796 0.476 0.6344
## Etnicidadefilipino -25.42861 19.91717 -1.277 0.2019
## Etnicidadehispanic other 4.25154 12.01266 0.354 0.7235
## Etnicidademexican 5.76263 4.63316 1.244 0.2138
## EtnicidadeOtherAsianOrPacific 0.72237 5.95321 0.121 0.9034
## Etnicidadeportoriq -8.50094 8.47990 -1.002 0.3163
## Etnicidadesouth americanIndian 12.94574 19.91416 0.650 0.5158
## Etnicidadewhite 3.57090 4.30742 0.829 0.4072
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.42 on 1389 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.2401, Adjusted R-squared: 0.2297
## F-statistic: 23.09 on 19 and 1389 DF, p-value: < 2.2e-16
## Analysis of Variance Table
##
## Model 1: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Marital + Gained +
## Smoke + Premie + Etnicidade
## Model 2: BirthWeightOz ~ Plural + Sex + MomAge + Gained + Smoke + Etnicidade
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1386 372384
## 2 1389 523802 -3 -151417 187.86 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
comparaCoeficientes(summary(fit_weight_6), sf1)Uma informação ausente sobre o peso do nascituro é: de qual bebê é o peso informado em BirthWeightOz? É a média dos dois ou três? É o menor (ou maior) deles? Essa falta de informação pode estar introduzindo alguma distorção no modelo. Façamos um sem gêmeos.
# fit_weight_7 é o fModelo 7 (sem gêmeos) ========
if (useEtn) { # Campo "etnicidade" unificado no modelo 7
fit_weight_7 = lm(
BirthWeightOz ~ Sex + MomAge + Gained + Smoke + Etnicidade,
data=ncb_weight)
sf7 = summary(fit_weight_7)
print(sf7)
anova(fit_weight_1, fit_weight_7, fit_weight_2) # compara variância de "1" com "2" e "7" =====
} else { # Campos étnicos separados (como no original) no modelo 7 =======
fit_weight_7 = lm(BirthWeightOz ~ Sex + MomAge + Gained + Smoke + RaceMom, data = ncb_weight)
sf7 = summary(fit_weight_7)
print(sf7)
fit_weight_7b = lm(BirthWeightOz ~ Sex + MomAge + Gained + Smoke + MomRace, data = ncb_weight)
sf7b = summary(fit_weight_7b)
print(sf7b)
fit_weight_7c = lm(BirthWeightOz ~ Sex + MomAge + Gained + Smoke + HispMom, data = ncb_weight)
sf7c = summary(fit_weight_7c)
print(sf7c)
anova(fit_weight_1, fit_weight_2, fit_weight_7, fit_weight_7b, fit_weight_7c)
}##
## Call:
## lm(formula = BirthWeightOz ~ Sex + MomAge + Gained + Smoke +
## Etnicidade, data = ncb_weight)
##
## Residuals:
## Min 1Q Median 3Q Max
## -99.731 -9.860 1.728 12.816 54.469
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 96.47907 5.47331 17.627 < 2e-16 ***
## SexF -2.13674 1.13286 -1.886 0.0595 .
## MomAge 0.45340 0.09778 4.637 3.87e-06 ***
## Gained 0.29870 0.04136 7.222 8.43e-13 ***
## SmokeY -7.56524 1.64242 -4.606 4.48e-06 ***
## Etnicidadeblack Black N -3.62693 4.78926 -0.757 0.4490
## Etnicidadeblack portoriq 6.08428 21.76065 0.280 0.7798
## Etnicidadeblack south american 4.24173 21.73921 0.195 0.8453
## Etnicidadecentro-south american 0.18372 6.49515 0.028 0.9774
## Etnicidadechinese 3.09240 15.72049 0.197 0.8441
## Etnicidadecuban 6.84812 15.71250 0.436 0.6630
## Etnicidadefilipino -23.97521 21.75047 -1.102 0.2705
## Etnicidadehispanic other 4.87315 13.11846 0.371 0.7103
## Etnicidademexican 4.93014 5.05936 0.974 0.3300
## EtnicidadeOtherAsianOrPacific 1.48282 6.50105 0.228 0.8196
## Etnicidadeportoriq -8.37477 9.26055 -0.904 0.3660
## Etnicidadesouth americanIndian 15.19172 21.74690 0.699 0.4849
## Etnicidadewhite 1.96807 4.70272 0.418 0.6756
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 21.21 on 1391 degrees of freedom
## (41 observations deleted due to missingness)
## Multiple R-squared: 0.0924, Adjusted R-squared: 0.08131
## F-statistic: 8.33 on 17 and 1391 DF, p-value: < 2.2e-16
## Analysis of Variance Table
##
## Model 1: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Marital + Gained +
## Smoke + Premie + Etnicidade
## Model 2: BirthWeightOz ~ Sex + MomAge + Gained + Smoke + Etnicidade
## Model 3: BirthWeightOz ~ Plural + Sex + MomAge + Weeks + Gained + Smoke
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 1386 372384
## 2 1391 625582 -5 -253197 188.48 < 2.2e-16 ***
## 3 1401 385454 -10 240128
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
comparaCoeficientes(summary(fit_weight_7), sf1) # , sf2, sf6)c. Verifique as premissas do modelo linear.
Os resíduos dos modelos 6 e 7 tem sua estatística W abaixo de 96%. Todos os outros ficara acima de 99%. Seus gráficos Q-Q apresentam deformidades que levantam suspeitas de que não sejam ~N.
# for (m in modelss) {
# qqnorm(m)
# }Peso do bebê
par(mfrow = c(2, 2))
plot(fit_weight_1)## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
par(mfrow = c(2, 2))
plot(fit_weight_1b)## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383
## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383
O gráfico de resíduos x valores ajustados apresenta uma forte heterodascidade, com uma grande concentração em volta de 123 Oz.
par(mfrow = c(2, 2))
plot(fit_weight_5)## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
Weeks
par(mfrow = c(2, 2))
plot(fit_weight_6)## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
Plural
par(mfrow = c(2, 2))
plot(fit_weight_7)## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
## Warning: not plotting observations with leverage one:
## 1200, 1306, 1383, 1387
d. Proposta de modelo
Com base nas análises, proponha um ou mais modelos lineares multivariados. Explique a sua escolha.
# bloco de código - item d# bloco de código - item e
weeks = c(10, # Poucas chances de sobrevivência
11,
23, # 17% de sobrevivência
25,
39, # termo
39.5,
40.2,
44,
45,
46, # risco
47)
gained = c(0.5, 11.3, 43.8, 100, 110)
age = c(4, 6, 8, 10, 12, 14, # Abaixo
21.5, 29.3, # no meio
44, 45, 46, 60) # acima
underage_premie = list(Plural="Single", MomAge=10, Sex="M", Weeks=20, Gained=5, Smoke="N")
predict(fit_weight_2, underage_premie)## 1
## 26.60699